AVERAGING ON MANIFOLDS BY EMBEDDING ALGORITHM 
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Abstract 

We will propose a new algorithm for finding critical points of cost functions defined on a dif- 
ferential manifold. We will lift the initial cost function to a manifold that can be embedded in a 
Riemannian manifold (Euclidean space) and will construct a vector field defined on the ambient 
space whose restriction to the embedded manifold is the gradient vector field of the lifted cost func- 
tion. The advantage of this method is that it allows us to do computations in Cartesian coordinates 
instead of using local coordinates and covariant derivatives on the initial manifold. We will exem- 
plify the algorithm in the case of SO{3) averaging problems and will rediscover a few well known 
results that appear in literature. 

MSC: 43A07, 53B21, 58A05. 

Keywords: Averaging, Distance functions. Rotations, Quaternions, Metriplectic dissipation. Gra- 
dient vector field. 

1 Introduction 

A problem which one frequently encounters in applications is to find an average for a finite set of sample 
points {yi, 2/2, • ■ • , J/r} belonging to a manifold N. A method to solve this problem is to construct a cost 
function Gn : A^ — > M associated with this set of sample points and then the average is defined by 

argminG7v(y)- 

Of special interest are the cost functions of least square type Gm{v) '■— X]I=i d.^{y:yi): where d is a 
distance function on N. The average of the set of sample points {j/i, 2/2, • ■ • , Vr} on the manifold N is 
the set defined by 



are min > d'' 

Z— 1 



iy,yi)- 



When the function Gjv is differentiable, this is equivalent with solving the equation dGN{y) = and 
with testing for which solutions the minimum value is attained. In order to write this equation we need 
the knowledge of a local system of coordinates on the manifold N, a requirement that might be difficult 
to fulfill in many practical cases. Another problem that we may encounter when this approach is that 
the set of sample points {j/i, ..., ?/,.} might not be entirely included in the domain of a single local system 
of coordinates. One way to overcome these difficulties is to lift the problem on a simpler manifold S. 

Let P : 5 — > A^ be a surjective submersion. The lifted cost function is Gs : S* -> E defined by 
Gs = Gn o P. The set equality P({x e S\dGsix) == 0}) = {y e N\dGN{y) = 0} shows that it is 
sufficient to solve the equation dGs{x) = on the simpler manifold S and project these solutions on 
the manifold N. 



A way to solve this new problem is to endow S with a Riemannian metric r and compute the critical 
points of the gradient vector field Vt-Gs and verify which one of them are local/global minima. We note 
that if xo G 5 is a critical point of the vector field WrGs, then it is a critical point of the vector field 
Vt'Gs, where r' is any other Riemannian metric on S. Although the manifold S might have a simpler 
geometrical structure then the initial manifold N, we still need a local system of coordinates or the 
knowledge of covariant derivatives on the manifold S in order to be able to compute the critical points 
of the function Gs, see [3], [3], [5], [TT], [TH], [H]. To further avoid the use of Riemannian geometry of 
the manifold S we will embed 5 in a larger space M. 

In the current paper we propose a solution for finding the critical points of the lifted cost function 
Gs by constructing a vector field vq £ X(M) on the ambient space M (usually an Euclidean space), 
that is tangent to the submanifold S an by also constructing a Riemannian metric t on S* such that 
Vols = ^tGs- Consequently, the critical points of the vector field vo (usually written in Cartesian 
coordinates) that belong to S are critical points of the lifted cost function Gs and their projection 
through the surjective submersion P gives the critical points of the initial cost function Gn- This 
construction will be illustrated in details for the averaging problem on the Lie group 5*0(3) associated 
with four different cost functions. 

In the literature the problem of averaging on Lie groups is often solved by using the exponential map 
as a tool to introduce local coordinates and so lifting the problem on the tangent space, see [3] , [1] , [H] , 
|12) . [13], [21]. This method can be generalized to the context of Riemannian manifolds as there also 
exists an exponential map. The exponential map can be further replaced by the more general notion 
of retraction developed in [3] , |10| . Averaging problems coming from real world applications have been 
studied in [2] , [5] , [H] , [23] , [H] , [55] by using the notions of covariant derivatives and the geometry of 
geodesies on various Riemannian manifolds. In this paper we will propose a different algorithm in order 
to deal with averaging problems on general differentiable manifolds. 

Other interesting cost functions are the ones that come from the Fermat-Torricelli problem. This 
averaging problem has been studied on various spaces, see [Z], [8], [15], [20,, and it will be of interest to 
apply the techniques developed in this paper to such problems. 

2 Averaging on a differentiable manifold 

We will solve the averaging problem presented in the Introduction by embedding the manifold S* as a 
submanifold of a Riemannian manifold {M,g). Further more, we will assume in this paper that S is the 
preimage of a regular value for a smooth function F := {Fi, . . . , Fk) : M — >■ M*"', i.e. S = F~^(co), for Cq 
a regular value of F. 

As we have already stated, our approach is to find a Riemannian metric r on the submanifold S and 
a vector field Vq € X(M) such that 

vo|s-V,Gs. (2.1) 

Consequently, our initial problem is equivalent with finding the critical points of the vector field Vq, 
which is defined on the ambient space M , and with choosing the ones that belong to S. 

The vector field that we are looking for is the standard control vector field vq introduced in [6] , which 
in turn is a continuation of the study we have begun in 0. Let {M,g) be a tti dimensional Riemannian 
manifold and Fi, . . . ,Fk,G : Af — > R be fc + 1 smooth functions. The standard control vector field is 
defined by 

vo = ^(-1)^+'=+^ det sj J;:;g^,,^^^^^)V^ + det sg;:S! VG, (2.2) 

4=1 

where det T,,^''"' p''. is the determinant of the Gramian matrix formed with the vector fields Vi^i, . . . , Vi^fc. 
The vector field Vq conserves the regular leaves of the map F :— (Fi, . . . ,Fk) : M — >■ M'^ and dissipates 
the function G, more precisely L^gG — detY^lp''"' ^'"^1 > 0. 

In what follows, we will describe the geometry of the standard control vector field. Let T : il^{M) x 



n\M) 



be the degenerate symmetric contravariant 2-tensor given by 

k 



5:(-ir+^-+Mets;j 






■.F,,...,F,) 
.,Fi,...,F^) 



VF, (E) VF, 



3 



deti];^::::£!g-. 



(2.3) 



The following result shows that vq looks like the gradient vector field of the function G, with respect 
to the degenerate symmetric contravariant 2-tensor T. 

Theorem 2.1. Jdjj On the manifold {M,g), the standard control vector field Vq is given by the following 
formula 

vo = idG(T). 

From the above theorem, computing the critical points of Vq (which is written in local coordinates 
on Af ) that belong also to S is equivalent with solving the following system oi m + k equations on M: 



F{x) = a 



Ti™(a;) 
T'"™(x) 



dG 
dG 

(x) 



dx'' 



(2.4) 



d 



d 



dG 



— — rdx'' . In general, a system of to + A: equations with to unknowns 



where T = TP«-— (g) -— and dG 

axP ax? 
might not have any solutions. But the above system that has m unknowns, namely the coordinates of 

X in M, has exactly to. functional independent equations due to the fact that the rank of the tensor T 

is TO, — A: for every regular point in the open set in M of regular points of F. 

The standard control vector field vo is tangent to every regular leaf and next we will recall its 

geometry when restricted to a regular leaf. Let ic '■ Lc —^ M be the canonical inclusion of the regular leaf 

Lc into the manifold M. The contravariant symmetric 2-tensor T is non-degenerate when restricted to 

a regular leaf Lc and consequently, the restriction can be inverted and it thus generates the Riemannian 

metric, see [6], Tc : ^{Lc) x X(Lc) -^ C°°(ic) defined by 



Theorem 2.2 

(i) Tc = 



detS 



(Fl 



re(X^ n := T-l((^,),X^ (z,),y-). 
6] On a regular leaf Lc we have the following characterizations. 






(ii) Vo|Lc = {ic)*^rAGoic). 

The Riemannian metric on S is t — Tcg and the above theorem also shows that the equality (12. ip 
holds, where Gs = G\s- 

In conclusion, in order to solve our initial averaging problem on the manifold N associated to the cost 
function Gjy, by using the standard control vector field, we apply the following embedding algorithm: 

(i) Choose the geometrical setting of a manifold S, a surjective submersion P : S* — > iV and construct 
the lifted cost function Gs ~ Gn ° P- 

(ii) Find a Riemannian ambient space {M,g) and a difFerentiable function F : M ^ M.'^ such that 
S = F~"'^(co), where cq is a regular value of F. Construct the contravariant symmetric 2-tensor T. 



(iii) Find a prolongation function G : M —> M. such that G\s 
vector field vq with the initial data F, G. 



Gs and construct the standard control 



(iv) Solve the system vo|5(a;) = (i.e. the system (|2.4p ) which is equivalent with finding the critical 
points of the lifted cost function Gs- Project through P these solutions and thus obtain the critical 
points of the initial cost function Gn- 



(v) Find critical points that are local/global minima for the cost fmiction Gn- 

Examples that fit the geometry underlying this embedding algorithm are given by the rich cases of 
quotient manifolds, N = S/K with K a group that acts on S. In the case when N can be written as a 
preimage of a regular value of a smooth map F, then (i) is not a necessary step of the algorithm {S — N 
and P is the identity). 

In the case when P*vo|5 is a vector field on the manifold N, then we have the set inclusion P({a; € 
5|vo|5(a;) = 0}) C {y G A^|P*vo|s(y) = 0}. If P is also a local diffeomorphism we will obtain 
equality between the above two sets and consequently, the critical points of the cost function Gn are 
characterized by the equation P*vo|5(y) = 0. In this particular setting the step (iv) can be replaced 
by the following: 

(iv') Suppose that P,vo|5 G X(A^) and P is a local diffeomorphism, then solve the equation P,Vo|5(y) = 
(the solutions are critical points for the cost function Gn)- 

For solving point (v) one can study the second order derivatives of the cost function Gn (see [S], 
|16) . |19) ) or study the convexity properties of the cost functions as in |15) . 

3 Averaging on S0{3) 

A problem frequently encountered in practice is finding a rotation that is in some sense an average of 
a finite set of rotations. The group of special rotations is not an Euclidean space, but a differential 
manifold. Consequently, an averaging problem on such a space has to be considered in the realm of 
differential geometry. 

The special orthogonal group 5*0(3) can be identified by a double covering map with the sphere S^ 
in M*, see [1], [23]. Our manifold S will be S^ and the ambient space (M, g) will be the Euclidean space 
R''. We will consider different distance functions on 5*0(3) that generate different cost functions and 
will solve the associated averaging problem, as we have described in the previous section, by using the 
embedding algorithm . 

From a historical perspective, we will use the quaternion notation on R^, although the quaternion 
structure will not be specifically used in our computations. The unit quaternions q = {q'^,q^,q'^,q'^) G 
S^ C K"' and — q E S^ C M^ correspond to the following rotation in 50(3): 

/ (g0)2 + (ql)2 _ (q2)2 _ (^3)2 2{q^q^-q°q^) 2{q^q^ + q°q^) 

R"! = 2{q^q^ + gOqS) (^0)2 _ (^1)2 ^ (^2)2 _ (^3)2 2{q^q^ - q^'q^) 

\ 2(glq3_^0g2) 2{q^q^ + q°q^) {q°? ^ {q^? - (q^ + (q^)^ 

This gives rise to the smooth double covering map P : 5^ — >■ 50(3), P(q) = R'". The covering map 
is a local diffeomorphism and consequently, instead of working with the distance function (cost function) 
on 50(3), we will work with cost functions on the unit sphere 5^ C M^ as it has been explained in the 
previous section. In this case, for the step (ii) of the embedding algorithm, we have S^ — F~^{1), with 
F : R" ^ R, i^(q) = (g°)2 + (gi)2 + (g2)2 _^ (^3)2^ rj^^^ tensor T is given by the formula, see [6] eq. 

(4.2), T = -VF(g) Vi^+ ||Vi^||2g-i, where g is the Euclidean metric on W^. Writing this explicitly, the 
matrix associated to the symmetric contravariant 2-tensor T is given by: 



[T(q)]=4 



-glqO (^0)2 ^ (g2)2 ^ (g3)2 

— g2q0 — g2gl (<Z°)^ 



q°q' 


-q"q' \ 


1 2 


1 3 


q q 


-q q 


q^f + {q'f 


-qS' 


q^q' 


{q'f + {q'Y + {q^Y ) 



In order to compute the standard control vector field vq according to the formula given in Theorem 
2.11 we need the following elementary computation: if aj(q) — ajo(q)(ig"+a;i(q)(i(7^+CL'2(q)rf'Z^+i^3(q)d'Z'^ 
is a one form on R^, then we have the vector field 

LT(q) = 4((q, q)cj(q) - (q, cj(q))q), (3.1) 



where we have made the notation aJ(q) := (a;o(q), Wl(q),|^2(q),'^3(q))• 
We will exemplify our construction by using various distance functions used in literature for mea- 
suring distances between Euclidean transformations. An extensive list is presented in [17 , where their 
equivalence and functional dependence is also studied. 

I. On the Lie group 5*0(3) we will consider the distance function di : 5*0(3) x 50(3) — ?> R+, 
di(Ri,R2) = ||Ri — R2||_F, where || • ||f is the Frobenius norm. The associated cost function is 
Gi,„,3, : 50(3) ^ R, 

r 

where Ri, ...,Rr are the sample rotation matrices, see [H], fFI\ . 

According to the step (i) of the embedding algorithm, we will lift the problem of finding the critical 
points of the cost function Gi^^.g, to the problem of finding the critical points of the lifted cost function 
Gig3 : 5^ -> R, Gig3 := Gig^^jj o P. 

In order to compute the lifted cost function Gi^g , by using the surjectivity of the covering map P, for 
each sample matrix R^ we will choose a sample quaternion q^ e 5"^ . We have the following computation: 

r r 

r 

= 2E(3-tr(R'i^R'i')) 

2=1 

r 
i=l 

where we have used the equality 

(q,q,)2 = i(tr(R'i^R'i-) + l). (3.2) 

For implementing the step (iii) we need to construct the prolongation function Gi : M^ — >■ K, 
Gi(q) = J2l=i 8(1 ~ (q^ qi)^)- Note that Gi^g and Gi are even functions so they do not depend on the 
choice of the sample quaternions q^ or — q^ which represent the same sample rotation R^. This subtler 
problem has also been addressed in a different way in pi^. Computing the differential 1-form dGi we 
obtain: 

dG^iq) = (^^(q),...,^(q)^ = (-IgJ^i" < ^^^^ >,■■■, -16 j^^' < q' ^^ >) 



j=l i=l 



= -16 E < q>q« > qj 



i=l 



Using (|3.ip . where the 1-form w is dGi, the system of equations (|2.4p becomes: 



(3.3) 



^^(q,q-*)((q,q)q-* - (q,q-*)q' 



j=i 

<q,q>=l 

and this represents the equation vo|53(q) — 0. The above system has four equations of third degree in 
the four unknowns q = (q*^ , q^ , q^ , q^) corresponding to the standard control vector field vq = plus 
the constraint equation that describes S^. 

The four equations corresponding to vq = are not functionally independent because, as we have 
stated in Section 2, the tensor T is degenerate. Nevertheless, since vo(q) G T^S^, for any q S 5^, the 
system p.3p can be described only by three equations plus the constraint equation. 



In order to write vo|s3(q) = as a system of three independent equations, we will push forward 
the vector field Vo|s3 trough the covering map P, which will also place us in the hypotheses of the step 
(iv') of the embedding algorithm. By definition of the push-forward operator, we have that 

P,vo|s3(R'») = TqP(vo|53(q)) e Tn.S0{3). 



As P is not an injcctive map, in order to obtain a vector field on the target space S'0(3), the equality 
T'qP(vo|s3(q)) = T_qP(vo|s3(— q)) must be satisfied. We will show in what follows that this is the 
case. 

We need to compute the tangent map of the covering map P : 5'^ — > 5*0(3). This map can be 
written as a restriction of the map P : M'^ -> R^ defined by P(q°, q^,q^, q^) = (((7°)^ + {q^f - (q^)^ - 
{q^)\ 2{q^q^ q"q^), 2{qV + q'q'), 2{q'q^ + q^q^), {q^f {q^f + {q^f {q')\2{q^q^ - q<^q^), 2{q^q^ - 



q^q^),2{q^q^ + q"q^), (g")"^ — (q^)^ — {q^Y + {q^^Y). where a rotation matrix R'* has been identified with 



a point in 



(the first line is identified with the first three components and so on) 



The matrix corresponding to the linear map TqP 



is given by 



JaCqP 



By direct computation, where we identify a vector in R^ with a tangent vector R'lA e T^q 5*0(3) 
{R'lA I A is a 3 x 3 skew symmetric matrix}, for q g S*"^ we have that 

TqP.Vo|s3(q)=TqP.Vo|s3(q) 

r 



= X!^^' q^) JaCqP • ((q, q)qj - (q, qj)q) 
1=1 

r 

^^(q,q,)R'lA,(q) 

r 

-R'l^(q,q,)A,(q), 



where 



A.(q) 



-,n„3 



-q^qf 



\„2 



T(il - q qt + 



q'q! 



q^q] - 
q^q" 



q\l 



-q'q! 



q'q] 



q^qf 



q'q"^ 



9 1 

q qt 



q^qf 



q'q'l 



q^qf 



-q'q] 



q'q!- 
^q'q^ 



9 n 

q q" - 

f q^qf 







q^ql 
- q^qf 



By identifying the vector JaCqP • ((q, q)qi — (q, qi)q) G K*^ with a 3 x 3 matrix, the skew-symmetric 
matrix is given by the formula 

A,(q) = Ri^JaCqP • ((q, q)q, - (q, q,)q). 



By noticing that Ai(— q) = — Ai(q), we obtain the equality TqP(vo|53(q)) = T_qP(vo[53(— q)) and 
consequently, P,vo|53 is a vector field on the manifold 50(3). 



As a result we obtain that the critical points of the cost function Gi 3 are the solutions of the 
following system (which is equivalent with 



^(q,q,)A,:(q)=0 _ ^^_^^ 

(q!q) = 1 

The advantage of the above system, compared with the system p.3l) , is that we have only three equations 
of second degree instead of four equations of third degree plus the constraint equation. 
If we transform the above system into rotations, by direct computation, we obtain: 

(q,q,)A,(q) = ^(R'l'^R'i - R'l^RiO (3.5) 

and the above system becomes: 

r 

^(Ri'^R*! Ri^R*!') ^ 0, 

2 = 1 

which is equivalent with 

R^Ri - Ri'^R = 0, (3.6) 

where R = 7 J2l=i ^^' = r ^l=i ^'- Solving this equation corresponds to step (iv') of the embedding 
algorithm and the solutions are the critical points of the cost function Gi^^.^. . Also, this equation is 
the same as the characterization for the Euclidean mean obtained in |19| . Thus, we have obtained that 
searching for critical points of the cost function Gig^.g, can be performed in two ways. More precisely, 
we can solve (|3.6p or we can transform the initial problem into quaternions and solve the system p.4[) 
and then transform the solutions into rotations. 

II. Next, we will consider the geodesic distance on the Lie group 5*0(3) which is defined by the 
angle of two rotations, see [T], g], [H], [H], gT], ^. For two rotations Ri,R2 G S0{3), d2(R-i,R.2) = 

||Log(Rf R2)||f — \/2\d\- where 9 £ (— 7r,7r) is the angle between the rotations Ri and R2. The 
associated cost function is 02^0^3) ■ SO{3) — > M, 

r 

G2,„<3,(R)=^||Log(RfR)||2,, 
1=1 

where Ri, ...,Rr are the sample rotation matrices. 

As before, we will compute the lifted cost function G2 3 as stated in step (i) of the embedding 
algorithm. More precisely, we havq^ 



92 



G2,3(q) - E l|Log(R^'^R'i)|||, = 2^ 

r 

-^ \ 2 



arccos 

1=1 



2 , tr(Ri'^Ri) - 1 



l t3T2t 

i=l 



2Earccos^(2(q,,q)2-l) 



2Earccos2(|(qi,q)|). 



i=l 



The lifted cost function G2 3 does not depend on the choice of sample quaternions q^ or — q^ that 
represent the same sample rotation R^. 



^\9i\ = arccos ( — '—= — ) , ii x a [—1, 1] then arccos(2a;^ — 1) = 2arccos(|x[). 



For step (iii) of the embedding algorithm, the prolongation function is given by G2 ■ M^\{0} 

21 I (q^ q« 



G2(q) — 2 y^ arccos^ 



|q||-l|q. 



The coefficients of the differential f-form dG2 are given by: 



dG2{ci)= (^(q),...,^(q) 



-4E 



sgn((q,q,))||q|pg^-|(q,q,)|g^ 



-4E 



rt l|q|IVIIqlPllqdP-(q>q.>' 

q,qj;i \ sgn((q,qj)) 



arccos 



i=l 



|q|| • m^ 



|q|IVIIqlPllqdP-(q'q.> 



l(q,q»)l 
|q||-l|q»l 

((q,q>qi: - (q,qi:>q)- 



In order to solve the differentiability problem we will eliminate from the domain of definition of the 
function G2 the hyperplanes 11^ = {q G M^ | (q, q^) = 0} and the lines di that go through points and 

In this case, the equation vo|53(q) = is equivalent with 

' V^ sgn((q,q,))arccos(|(q,qj)|) \ ^ n 

2^ ,, , ,, ((q, q)q» - (q, q»)q) = o 

r 

q e ^MJ n, 

I qT^iqi 

By using the same arguments as before, the equation vo[53(q) = is equivalent with P»vo[53 — 
which has the following form 



y^ sgn((q,q,))arccos(|(q,q,)|) ^ 
^ Vl-(q,q.>2 

r 

q e ^MJ n,; 

i=l 

I q 7^ ±qi 

The above expression also shows that P*vo|53 is a vector field on the manifold 5*0(3). 
In order to transform (|3.7p into rotations we need the following computation: 



(3.7) 



sgn((q, q.,) ) arccos(|(q, q^) |) arccos(|(q, q,) |)(q, q., 



d 



Vi-(q,q^)' 



l(q,q»)IVi-(q>q»)' |cos(f )|y^i-cos2(|) 



iq,q» 



sm 



(q, q* 



sin(0O 



(q,q»;, 



where we have used the property q £ S^\ lJi=i ^i which implies that (q, q^) 7^ 0. Using p.5p we obtain 
the equivalent system in rotations, and this corresponds to the step (iv') of the embedding algorithm. 



V(Ri'^Ri - Rq^R^')— -i- = ^ VLog(Ri'^R'i) = 0. 



sm 



(3.; 



The above equation has been also obtained in 19 as a characterization of the Riemannian mean. This 
was to be expected, as we are averaging the same cost function by two different methods that naturally 
lead to the same result. Also note that the cost function G2S013) i^ not well defined for angles 9i = ±7r 
or equivalently (q, qi) = 0, a situation in which the prolongation function G2 is not differentiable. 



III. Another distance used on the Lie group S0{3) is given by d3(Ri,R2) = 1 — ^•\/tr(RfR2) + 1. 
This distance appears in |17) , where the functional dependence between d^ and d2 is also given. We 
notice that this distance function on 5*0(3) can be obtained from the pseudodistance da : S*^ x S*"^ — > M, 
(i3(qi, q2) = 1 — |(qi, q2)|- The cost function Gs^^.^. (R) = X]i=i '^K^^i ^^0 i^ lifted to the cost function 

r 

The prolongation function is G3 : M'' — > M, G3(q) = X]I=i(l ^ I(QjIj)I)^- The coefficients of the 
differential 1-form dG^ are given by: 

r 

dG;(q) =-2^(1- |(q,q,)|)sgn((q,q,))q„ 

i=l 

where, in order to have differentiability, we restrict our selves to the open set q e S^\ lJl=i ^i (the set 
Hi is the hyperplane determined by qi as before). 
The system of equations (|2.4p becomes 

r 

X!(l ^ l(q,qz)l)sgn((q,q»))((q,q)qj - (q,qj)q) == 

*=^ . (3.9) 

qe Aljn^ 

and this represents the equation Vq = restricted to the open set S^\ lJl=i n^ of the sphere S'^. 

As in the previous cases, after applying the projection operator P*, we generate the vector field 
P*vo|53 and consequently, we obtain the equivalent system of equations: 

jr r 

^(l-|(q,q,)|)sgn((q,q,))A,(q)=0 

'=1 r ■ (3.10) 

qe53\Un, 

1=1 

Transforming the above system into rotations, after some algebraic manipulations using (13. 5p and 
, we obtain the equation from the step (iv') of the embedding algorithm 

J2 ^ 1 (R^'^R'' - R^'^R''-) = 0. 

,;=i \ ytr(Rq^RqO + 1 / 

IV. Next we will apply the techniques developed so far for a cost function that is of L^-type but 
not of least square type, see [HI- We start with the distance di : 5*0(3) x 50(3) -^ M+, (ii(Ri, R2) = 
||Ri — R2||_F. The associated L'*-mean cost function is 043^,3, : 50(3) -^ M, 



^4so(3) (R) - 2J I 1^ ^ ■^■ 



4 
'i||_F- 



1=1 



The lifted cost function is given by 64^3 : S^ 



\2\2 



G4,3(q) = 64^(1 -(q,q- 

1=1 

The prolongation function is G4 : R"* ^ M, G4(q) == 64X;l=i(l - (q,qi)^)^ 



As in the case of the function Gi, the system (|2.4p corresponding to the function G4 becomes: 

r 

Y^{{i,<\i) -(q,q»)^)((q,q>qi-(q,q»)q) = o 

<q,q>=l 

By using the same arguments as before, in the case of function G4, the equation Vo|53(q) = is 
equivalent with P^VQigs = 0, which has the following form 

r 

^(l-(q,q,)2)(q,q,)A,(q)=0 ^ ^^ ^2) 

(qjq) = 1 

Using again (I3.2p and (|3.5I) we obtain, according to step (iv') of the embedding algorithm, the equation 
in rotations that gives the critical points of the cost function Gi^^.^. , 



^(3 - tr(Ri'^RiO)(R-''''^R'' - R^^R"') = 0. 



Acknov^rledgements. Petre Birtea and Dan Comanescu have been supported by the grant of the 
Romanian National Authority for Scientific Research, CNCS UEFISCDI, project number PN-II-RU- 
TE-2011-3-0006. 

References 



[4: 



[9 



[10 



S.L. Altman, Rotations, quaternions, and double groups, Oxford University Press, 1986. 

P.A. Absil, R. Mahony, R. Sepulchre, Riemannian geometry of Grassmann manifolds with 
a view on algorithmic computation, Acta Applicandae Mathematicae, Vol. 80, Issue 2 (2004), pp. 
199-220. 

P.A. Absil, R. Mahony, R. Sepulchre, Optimization Algorithms on Matrix Manifolds, Prince- 
ton University Press, 2008. 

V. Arsigny, X. Pennec, and N. Ayache, Bi-invariant means in Lie groups. Application to 
left-invariant polyaffine transformations, INRIA research report, 5885 (2006). 

P. Birtea, D. Comanescu, Asymptotic Stability of Dissipated Hamilton-Poisson Systems, SIAM 
J. Appl. Dyn. Syst., Vol. 8, Issue 3 (2009), pp. 967-976. 

P. Birtea, D. Comanescu, Geometric Dissipation for dynamical systems, Comm. Math. Phys., 
Vol. 316, Issue 2 (2012), pp. 375-394. 

D. Comanescu, S.S. Dragomir, On the Torricellian Point in Inner Product Spaces, Demonstra- 
tio Mathcmatica, Vol. XLI, Issue 3 (2008), pp. 639-650. 

D. Comanescu, S.S. Dragomir, E. Kikianty, On the Torricellian Points in Normed Linear 
Spaces and Applications, (2012), http://ajmaa.org/RGMIA/papers/vl5/vl5a3.pdfi 

A. Edelman, T.A. Arias, S.T. Smith, The geometry of algorithms with orthogonality con- 
straints, SIAM J. Matrix Anal. Appl., Vol. 20, Issue 2 (1998), pp. 303-353. 

S. Fiori, Lie- group-type neural system learning by manifold retractions. Neural Networks, Vol. 21 
(2008), pp. 1524-1529. 



10 



[11] S. Fiori, Solving Minimal- Distance Problems over the Manifold of Real-Symplectic Matrices, SIAM. 
J. Matrix Anal. & AppL, Vol. 32, Issue 3 (2011), pp. 938-968. 

[12] S. Fiori, T. Tanaka, An algorithm, to compute averages on matrix Lie groups, IEEE Trans. Signal 
Process., Vol. 57, no. 12 (2009), pp. 4734-4743. 

[13] V. M. Govindu, Lie-Algebraic Averaging For Globally Consistent Motion Estimation, Proe. IEEE 
Conf. on Computer Vision and Pattern Recognition, Washington, DC, Vol. 1 (2004), pp. 684-691. 

[14] C. Gramkow, On Averaging Rotations, J. Math. Imaging Vis., Vol. 15 (2001), pp. 7-16. 

[15] R. Hartley, J. Trumpf, Y. Dai, Rotation averaging and weak convexity. International Sympo- 
sium on Mathematical Theory of Networks and Systems, Piscataway, USA (2010), pp. 2435-2442. 

[16] U. Helmke, M.A. Shayman, Critical Points of Matrix Least Squeres Distance Functions, Linear 
Algebra and its Applications, Vol. 215 (1995), pp. 1-19. 

[17] D.Q. Huynh, Metrics for 3D Rotations: Comparison and Analysis, J. Math. Imaging Vis., Vol. 
35 (2009), pp. 155-164. 

[18] F.L. Markley, Y. Cheng, J.L. Crassidis, Y. Oshman, Averaging Quaternions, Journal of 
Guidance, Control, and Dynamics, Vol. 30, Issue 4 (2007), pp. 1193-1197. 

[19] M. Moakher, Means and averaging in the group of rotations, SIAM J. Matrix Anal. AppL, Vol. 
24, Issue 1 (2002), pp. 1-16. 

[20] P.L. Papini, J. Puerto, Averaging the k largest distances among n: k-centra in Banach spaces, 
J. Math. Anal. AppL, Vol. 291, Issue 2 (2004), pp. 477-487. 

[21] X. Pennec, Computing the mean of geometric features - Application to the mean rotation, Institut 
National de Recherche en Informatique et en Automatique, Rapport des Recherche 3371, 1998. 

[22] I. Sharf, A. Wolf, M.B. Rubin, Arithmetic and geometric solutions for average rigid-body 
rotation. Mechanism and Machine Theory, Vol. 45 (2010), pp. 1239-1251. 

[23] M.D. Shuster, A Survey of Attitude Representations, Journal of the Astronautical Sciences, Vol. 
41, Issue 4 (1993), pp. 439-517. 

[24] C. Samir, P.A. Absil, A. Srivastava, E. Klassen, A Cradient-Descent Method for Curve 
Fitting on Riemannian Manifolds, Foundations of Computational Mathematics, Vol. 12, Issue 1 
(2012), pp. 49-73. 

[25] R. Subbarao, P. Meer, Nonlinear Mean Shift over Riemannian Manifolds, Int. J. Comput. Vis., 
Vol. 84 (2009), pp. 1-20. 

[26] B. Vandereycken, P.A. Absil, S. Vandewalle, A Riemannian geometry with complete geodesies 
for the set of positive semidefinite matrices of fixed rank, IMA Journal of Numerical Analysis, (2012), 
doi:10.1093/imanum/drnxxx 



11 



